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Abstract 

We develop a general form of the Ritz method for trial functions that do 
not satisfy the essential boundary conditions. The idea is to treat the latter 
as variational constraints and remove them using the Lagrange multipliers. In 
multidimensional problems in addition to the trial functions boundary weight 
functions also have to be selected to approximate the boundary conditions. We 
prove convergence of the method and discuss its limitations and implementa- 
tion issues. In particular, we discuss the required regularity of the variational 
functional, the completeness of systems of the trial functions, and conditions 
for consistency of the equations for the trial solutions. The discussion is accom- 
panied by a detailed examination of examples, both analytic and numerical, 
to illustrate the method. 

Keywords: Convex functional, boundary value problem, essential bound- 
ary condition, removing variational constraints, energy space, minimizing se- 
quence, trial function, complete system, convergence of trial solutions 
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1 Introduction 



In variational problems linear boundary conditions are often divided into essential 
(geometric) and natural (dynamic) [H 11.12], [21 4.4.7]. More generally, one calls the 
boundary conditions essential if they involve derivatives of order less than half of the 
order of the differential equation, and natural otherwise [3j 1.1.2]. The common lore 
on the Ritz method is that the trial functions may violate the natural conditions, but 
must satisfy all the essential ones [2] 4.4.7], [4]. The reason is that the variational 
equations force the natural conditions on the trial solutions anyway, even if the trial 
functions themselves do not satisfy them. 

But what if we wish to use trial functions that violate the essential conditions 
as well? For instance, in problems involving parametric asymptotics the trial func- 
tions are pre- imposed with no regard for boundary conditions [HI E] , and in initial- 
boundary problems with time-dependent boundary conditions the (time indepen- 
dent) trial functions can not satisfy them in principle. One may also wish to use 
such violating trial functions because they are simpler. Sure, it is easy enough to 
adjust them in one-dimensional examples, but it is not so easy at all in higher di- 
mensions, especially in problems with fancy boundaries. Thus, there is abundant 
motivation to generalize the Ritz method to the trial functions that do not satisfy 
the essential conditions. However, surprisingly few authors consider such generaliza- 
tions, many works on the subject are rather old, and sometimes give prescriptions 
that produce inconsistent systems and/or erroneous approximations [7J [8]. 

This is not to say that nothing has been done at all. One idea is to treat the 
essential boundary conditions as variational constraints and to remove them as any 
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other constraints using the Lagrange multipliers. This idea is so natural that it 
appears occasionally in some applied works at least since 1946, see [9], where the 
authors explicitly cite the simplicity of the trial functions that violate the essential 
conditions as a reason for using them. Violating trial functions were even used in the 
boundary eigenvalue problems for vibrating plates, see [10] and references therein. 
However, in all of these works Lagrange multipliers are essentially applied by analogy 
to finite dimensional optimization, numerical success serving as the only validation. 
A systematic discussion, or even a general description of the method, accounting for 
its requirements and limitations seems to be missing in the literature. This is all 
the more surprising since, as we shall see, the standard theory of the Ritz method 
[TT1 [T2l IT3] can be readily accommodated to handle such a generalization. 

A very interesting paper [13] discusses one particular issue that arises when using 
the Lagrange multipliers in one-dimensional boundary value problems, completeness 
of the trial functions. We will discuss this and other issues that come up, extend 
the method to higher dimensions and give a general proof of its convergence. Our 
analysis indicates that more care is needed compared to the usual Ritz method, 
to wit, the variational functional has to be more regular on a larger space, the trial 
functions have to be complete in this larger space as well, and the multipliers can not 
be eliminated from the approximating systems using the usual variational formulas 
because of convergence issues. In the higher dimensional problems one needs to 
select boundary weight functions in addition to the trial functions, and balance the 
numbers of each to obtain well-posed approximating problems. 

We try to keep our discussion more suggestive than technical, so the paper is 
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structured somewhat unusually. We start in Section [2] by presenting some simple 
one- dimensional examples, and apply to them some natural looking approximating 
procedures recommended in the literature when the trial functions do not satisfy the 
essential boundary conditions. The examples are so simple that not only the exact 
solution, but even all trial solutions can be computed analytically. Distressingly, in 
many cases the approximations do not exist, do not converge, or converge to a wrong 
answer. The reasons turn out to be subtle, but they all can be traced to the vagueness 
in the concept of " approximation" . In Sections El H] we flush out the false assumptions 
behind the failed approximations and develop the Ritz-Lagrange method that avoids 
them, proving along the way that it works. Unfortunately, in its original form the 
method only applies to one-dimensional problems, so more developing and proving is 
done in Section [5] for higher dimensions. Numerical applications to multidimensional 
problems follow in Section El The paper ends with Conclusions, where we summarize 
our findings and discuss Galerkin type generalizations. Technical proofs are collected 
in the Appendix. 



2 One-dimensional examples 

In this section we test several reasonably looking methods on simple one-dimensional 
problems. Not all of them work and, as we will find out later, some of them work 
for the wrong reasons. These examples will highlight some subtleties and serve as 
a motivation for developing a general method later on. A typical procedure for 
numerical solution of boundary value problems is to represent the trial solution as 
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a linear combination of finitely many trial functions and solve a finite dimensional 
problem that results. Most frequently the trial functions are required to satisfy the 
essential boundary conditions, otherwise some additional effort is needed to take care 
of them. 

Problem 1. Consider a boundary value problem for the second-order equation 
u X x = f on [0,L] with essential boundary conditions on both ends of the interval 
u(0) = u(L) = 0. We set for convenience L = n, and / = 1 to make everything 

1 7T 

explicitly computable. The exact solution is easily found to be u = -x 2 — —x. 

Lanczos tau method. This is perhaps the most popular method that uses 
trial functions not satisfying essential conditions. The reader is warned that we use 
the broad understanding of the method, as in |5] [6] (review [7] also describes it 
without using the name). Many authors, including Lanczos himself, understood it 
much more narrowly, only allowing Chebyshev polynomials as trial functions, see 
[T5] . The basic idea is to make the residual, the difference between the left and the 
right hand sides of the equation, orthogonal only to some trial functions used in 
the trial solution. This creates a shortfall in the number of equations compared to 
the number of unknown coefficients. This shortfall is used to impose the essential 
boundary conditions directly. 

In our case the system to be solved can be written succinctly as 



where selected trial functions are substituted for the weight Su. We select cos nx, n > 




(1) 
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as our trial functions, they obviously do not satisfy the boundary conditions. Tak- 

N 

ing N of them the trial solution is of the form = h a n cos nx with un- 

2 n=l 

known coefficients Oj (| in front of is for agreement with the convention for the 
cosine series). Since there are two boundary equations we need to omit two trial 
functions when choosing weights, e.g. take 8u = cosmx, m = 0, . . . , N — 2 in Eq.([T]). 
This gives 

/ (u xx — l)Su dx = = / (— n 2 a n cos nx — 1 j cosmx dx = 0, (2) 
Jo Jo v n=1 J 

where m = 0, . . . , N — 2. However, these equations are already inconsistent for any 
N since J Q cosnxdx — for n > 1 and the m — equation becomes J — 1 dx — 0. 
If we choose a different subset of trial functions as weights, say m — 2, . . . , N, then 
the equations are consistent, but we find — m 2 a m | = 0, i.e. a m = for m > 2. The 
boundary conditions then imply that also a® = ai = yielding = for all N. 
One can see that any choice of weights here leads to inconsistency or to the trivial 
solution. The Lanczos tau method fails completely for our choice of trial functions. 

Ritz method with boundary terms. Let us turn to the Ritz method now. The 
corresponding variational functional is J(u) = \ {u x ) 2 + fudx and the boundary 
value problem is equivalent to minimizing it on functions satisfying the boundary 
conditions. Since our boundary conditions are essential, and our trial functions do 
not satisfy them, the standard approach would have to be modified in some way. It 
turns out that not every plausible modification works. One obvious idea is to treat 
the essential conditions as variational constraints and remove them using Lagrange 
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multipliers. The Lagrange functional is £ = J + Xqu(0) + A 7r u(7r), where Ao, A„- are 
the Lagrange multipliers, and the variational equation is 

= 52, = 5 J + \ Q 5u(0) + \ n 5u(7r) + u(0)5X + w(7r)5A 7r 

u x Su x + fSudx + X 5u(0) + X n Su(ir) + u(0)5X + M(7r)5A 7r 

(— u xx + f)8u dx + u x 5u + Xo5u(0) + \ n 5u(7r) + u(0)5\o + ■u(7r)5A 7r 

o 

(3) 

(—u xx + f)5u dx + ^A,,- + u x (n) j5u(tt) + ^Ao — u x (0)j5u(0) + u(0)5\o + u(ir)5\ 7T 

Since the boundary variations are independent of the internal ones we see that Ao = 
u x (0) and X n = —u x (tt). This allows us to eliminate the multipliers from Eq.Q 
leading to: 

/'7T 

0= / (-u xx + f)5udx + u(0)5u x (0)-u(7r)5u x (TT). (4) 



JO 

This equation is derived as a generalization of the usual Ritz system e.g. in [8]. Let 
us test it on our example. As before, one substitutes the trial solution for u 
and the trial functions for the weight 5u. With our trial functions cos mi however 
5u x (0) = 5u x (tt) = for all m. But then Eq.fllJ) reduces to the first equation in 
Eq.Q and produces the same system as in Eq.(j2J), only with m = 0, . . . , N. We 
already saw that the m = equation can not be satisfied, making the entire system 
inconsistent for any N. This method does not work either. 

Ritz-Lagrange method. Perhaps, instead of trying to eliminate the Lagrange 
multipliers from Eq.flBJ) we should keep them as unknowns and supplement the ob- 



8 



tained equations with the boundary conditions, just as one would do for any vari- 
ational constraints. One positive trait is that the number of unknowns will always 
match the number of equations since each multiplier corresponds to a boundary 
condition. Instead of relying on Eq.Q let us find £(m^) explicitly. Recall that 

N 

(N) a . 

ir ' = h a n cos nx, so 



n=l 



N / N \ / JV 



£ (>)) = |a + \ £ -X + Ac | + £ a J + A. f + ^(-l)"a n ) . (5) 

n=l \ n=l / \ n=l 

Variational equations are 7: — = and adding the boundary conditions we get the 
following system: 



^0 \ — v 

7r + A + A 7r = — + 2_^a n = 

n=l 

2 N 

^-a n + A + (-1)^ = ^ + ^(-l) n a n = . 



There are N + 3 equations for iV + 3 unknowns ao, • • • , ojv, Ao and A^. For even n^O 
one immediately finds that a n = \ from the first column equations. Analogously, 
for odd n we have a n = — ^(Xo ~ A^). Then subtracting the second equation in the 
second column from the first 

J2 2a n = --(\ -X 7T ) V^ = 0. 

n odd n odd 

This implies that \o = X n = an d all the odd coefficients vanish. Thus, we 
have a n = — ni for n > 1 and only ao remains undetermined. Adding the second 
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column equations 

a + 2a n = a + — - = , and therefore 

n even n even 

LA72J iN/2] 2 

a ° = ^ >:= - 4 £(^ = -Ep ^ -y ' 

where ['J is the floor function returning the largest integer not exceeding its argu- 
ment. Summarizing, we conclude that the trial solutions converge to the sum 
of the series 

7T 2 °° 1 _|_ / 1 V" 

u Wf x ) yu^\x) :=-— + V \ J cosnx. 

v 7 jv-^oo v ; 12 ^ n 2 

n=l 

Recall that by extending a square integrable function w(x) on [0, ir] to an even func- 

2 Z" 71 

tion on [— tt, it) one can expand it into a cosine series with coefficients a n = — / it;(x) cos nx dx 
[T6t 12.1]. Considering that 



7T 2^1-(-l) n . 2 7T 2 ^(-1)™ 

> 5 cosnx and x = h 4 > — cosnx (7) 

2 7r ^ n 2 3 ^ n 2 

n=l n=l 



we see that u^ 00 " 1 is exactly the cosine series for the exact solution u = -x 2 x. 

3 2 2 

Finally, we have a method that works for this example at least. To make sure it 
was not an accident we will now apply it to some other boundary value problems. 

Problem 2. Consider the biharmonic equation u xxxx = f with the boundary 
conditions -u(O) = u(ir) = u xx (0) = u xx (it) = 0. For / = 1 the exact solution is 

4 s 
U = 7T h 7T . 

24 12 24 
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We apply the Ritz-Lagrange method again. The last two boundary conditions 
are natural and we need not worry about them. Thus, the variational problem 



is to minimize J = f Q 



L 1 



fu dx on the space of functions satisfying the 



essential boundary conditions only. Now let us find the trial solutions. Taking 

N 



U 



(TV) _ ^0 



+ a n cos nx we have 



n=l 



rfu (fl) ) = J + \ u(0) + A^(tt) = 

N f N \ ( N \ 

n=l \ n=l J \ n=l J 



TV 7T 

~2 a ° 



since cos nx = (— l) n . The Ritz-Lagrange system is 



— 7r + A + A„- = 



N 



irn- 



a n + X + {-l) n K = ^ + ^T(-l) n a„ = 



71=1 

N 



(9) 



n=l 



It can be solved along the same lines as system Eq.((6]) and we get Ao = A^ = |, 
a n = - 1+( ~4 1)U for n > 1 a = aj^ := J Y}^J? jk — > where [-J is the floor 



N->oo 



function returning the largest integer less than or equal to its argument. Using Eq.flTJ) 
and 



7T 



x 



+E( 67r 



n=l 



■1^ 



12 1 



n 



n rr 



7T 



cos nx ; x 



n=l 



J-l) n ( 

?tt 2 ^ — 48- 



rr 



rr 



(10) 
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one can check that 

n=l 

and therefore w — li*- 00 -* = n 3 tt 2 — . This time the trial solutions do not converge 

24 24 & 

to the exact one, and the limit difference is a linear combination of x and x 2 . 

Problem 3. Consider the biharmonic equation u xxxx = f again, but now with 
all essential boundary conditions u(0) = u(ir) = u x (0) = u x {ix) = 0. For / = 1 the 

4 3 2 

CC tJC kJO 

exact solution is u(x) = tt h tt 2 — . 

v ; 24 12 24 

The Lagrange functional will now be £ = J + A w(0) + A 7r w(7r) + \' u x (0) + A^-u(7r) , 
but if we use our trial functions cos nx, n > the last two terms in £ will be for any 
v,( N '. But then £ is the same as in Eq.([H]), and therefore the Ritz-Lagrange system is 
the same as in Eq. ([9]) . Then the trial solutions are also the same and they converge 
to from Eq- ffTT]) . However here, unlike in Problem 2, this is the right solution. 
The reader may wish to entertain herself by thinking over the last two problems. To 
help we will remove one red herring, as we shall see the presence of natural conditions 
is not an issue here. 

3 Continuity and convergence 

In this section we will analyze what went wrong (and right) in our examples and 
come up with a general scheme that provably works. The Lanczos tau method is the 
easiest to figure out. The idea behind the method is fairly intuitive. Since the exact 
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solution has residual it is certainly orthogonal to all trial functions, and it is the 
only such function that also satisfies the boundary conditions. So as we construct 
approximations with residuals orthogonal to more and more trial functions, while 
also forcing the boundary conditions on them, it stands to reason that they should 
approach the exact solution in some sense. There are several assumptions that this 
argument relies on however. First, the system of trial functions should be complete, 
i.e. one should be able to approximate any function by their linear combinations. 
Second, we have to make sure that in the limit the residuals do become orthogonal 
to all trial functions. This was not the case in our second application of the Lanczos 
tau method, when cosmi, m = 2, . . . , N were used as weights. Indeed, no residual 
ever had to be orthogonal to 1 or cosx. This is how the spurious solution u = 
slipped through the cracks with the residual —1. 

It is harder to explain why we got inconsistent systems when using cos mx with 
m — 0, . . . , N — 2. Let us use the benefit of hindsight and look at the cosine series 

yj-2 ~ I _|_ (-1)™ 

for the exact solution u — h > t, cos nx. Cutting it off at N terms 

12 n 2 

n=l 

and computing the residual we get u xx — 1 = — Y^n=\ (l + ( — -0") s ^ nnx - Clearly, 
this series does not converge to in any apparent sense. But our argument above 
relied on exactly this kind of continuity assumption: if trial solutions approach the 
exact solution so do the residuals, and this is not the case here. Without the exact 
solution to approach 'trial solutions' have nothing to approximate, hence there is no 
reason for the Lanczos tau system to be consistent. And indeed it was not. 

Looking closely at the above reasoning the reader will notice the vagueness in 
the meaning of "approximation". Approximation in what sense? To talk about ap- 
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proximating we need a way to measure distance between functions. Such measure 
is usually provided by Banach space norms. The norm relevant to the discussion 
above is the L 2 norm \\u\\l 2 := {^J \u\ 2 5udx^j with the corresponding inner prod- 

uct (if, v) := / uv dx. The Banach space of square integrable functions with the L 2 
Jo 

norm is a Hilbert space is denoted L 2 ([0, 7r]), and all above references to convergence, 
completeness, continuity and orthogonality referred to this space. 

But the L 2 space is not sufficient for a discussion of variational methods. Many 
variational functionals come with a natural space of their own, the energy space [1, 
Ch.2]. For the functional of Problem 1 the energy space is W^QO, %]), the Hilbert 
space of functions square integrable with their first derivatives and vanishing at 
and 7r, with the norm : = (|I m I|l 2 + II^IlL) 2 • This norm ' s stronger than 

the L 2 norm in the sense that any W 2 convergent sequence converges in L 2 , but not 
conversely. On its energy space a variational functional J typically has two important 

properties: it is continuous and it is growing at infinity, i.e. J{u) > 00 [121 

||m||— >oo 

III. 10. 2]. For a convex functional (we will only consider those) these two properties 
are sufficient to prove convergence of the usual Ritz approximations in the energy 
norm [HI 6.2A]. Note that continuity and growth conditions balance each other: for 
a stronger norm continuity is preserved of course, but the functional may no longer 
grow at infinity, and vice versa for a weaker norm. 

For our purposes the concept of energy space is not quite suitable because it 
hardwires essential boundary conditions are hardwired into it, and our trial functions 
do not satisfy them. Instead we start with a reflexive Banach space U (the reader will 
not lose much by assuming it to be Hilbert), that has nothing to do with the boundary 
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conditions, and a convex functional on it J : U — > M. Next, we introduce the 
boundary operator, a linear map r : U — > W that maps functions into their boundary 
values. The subspace U := {u £ U \Vu = 0} consists of functions that satisfy the 
boundary conditions. For Problem 1 we have U = W^QO,^]), Tu = (w(0), u(7r)) T 
and 1A = W^dP) 7r]). The following three assumptions turn 1A into a generalized 
analog of the energy space: 

1. J : U — > R is convex and continuous; 

2. J(u) > oo, i.e. J grows at infinity, on U; 

\\u\\— >co 

3. r : W — >■ M s is linear and continuous. 

This setup applies to homogeneous boundary conditions only. Non-homogeneous 
boundary conditions can be accommodated in the usual manner, by selecting a func- 
tion that satisfies them and switching to the differences with it. They solve the 
corresponding homogeneous problem and all convergence issues can be reduced to 
it, see e.g. [17, 2.1]. 

The functional and the boundary conditions being dealt with we now turn to 
the trial functions. Recall that a system of elements in a Banach space is called 
complete if any element has their linear combination within any given distance from 
it. Let {(pi} be a complete system in hi and let denote the linear span of 

cj>i, . . . , 4>n- The Ritz-Lagrange approach to approximating the minimizer of J on U 
amounts to minimizing it on IA^ N > subject to the boundary conditions. This is of 
course equivalent to minimizing it on := U^ N ^ nU, whose elements are such 

linear combinations of the first iV trial functions that satisfy the boundary conditions. 
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Although by assumption about completeness all functions in U can be approximated 
by linear combinations of the trial functions, it is not a priori clear that functions 
from U can be approximated by linear combinations that are themselves in U. The 
next lemma proved in the Appendix assures us that this is the case. 

Lemma 1. For any complete system of elements in U there exists a system of their 
finite linear combinations belonging to U. which is complete in U. 

This lemma effectively reduces the Ritz-Lagrange method to the traditional Ritz 
method. Indeed, if {fa} is the complete system in U produced by Lemma [T] then 
applying the Ritz method with (f>i as the trial functions instead of fa amounts to 
minimizing JonW ( N >. In other words, the Ritz-Lagrange method with {fa} produces 
the same (up to re-indexing) as the Ritz method with {4>i}. This allows us to 
use well-known results on convergence of the Ritz method [T2l IV. 12], (TTJ 6.2A] to 
prove convergence of its Ritz-Lagrange generalization. 

One additional issue left is that one needs J to be differentiate to minimize it 
using Lagrange multipliers. Note that we only need this differentiability on finite 
dimensional subspaces U^ N \ that was explicitly true of J in all our examples. In 
general, it suffices that J is Gateaux differentia ble on U, i.e. that its first variation is 
linear in 5u p2} 1.2.1], [TT1 3.2]. The differentiability also automatically guarantees 
continuity. For general convex functionals only weak convergence can be expected, 
we use symbol — > to denote it. However, due to Sobolev embedding theorems [TTJ, 
1.6], [18,, 1.8] weak convergence in W\ for example implies convergence by norm in 
L 2 . Summarizing we get the theorem below. The proof is fairly standard, but we 
outline it in the Appendix for the convenience of the reader. 
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Theorem 1. Suppose J :IA — )■ R is Gateaux differentiable, T : U — > W is a bounded 
linear operator, and {4>{\ is a complete system in U. If J is convex and grows at 
infinity on U then it has a minimizeru on it, as well as minimizers on allU^ N \ 
and there exists a subsequence iV& such that Moreover, if J is strictly 

fc— >oo 

convex on U then u, u ( N > are unique and u ^ — — y u. In both cases the values of J 

fc— ^oo 

converge to its minimum on U. 

One can say more for quadratic functional of the form J(u) = ^B(u,u) + l(u), 

where B is a symmetric bilinear form and / is a linear form. These types of 

functionals produce linear boundary value problems [HI 22.1]. In Problem 1 we 

had B(u,v) = / u x v x dx and l{u) = / fudx, while in Problems 2, 3 they were 
Jo Jo 

B(u,v) = / u xx v xx dx and l(u) — — I fudx. Such J are convex if B is positive 



definite, and strictly convex if B is strictly positive definite, i.e. B(u,u) > £||u|| 2 for 
some e > 0. One can check that our functionals are strictly convex on U in both 
cases [TFJ 1.8]. Simple algebra shows that for any u,v GW 

J{u) - J{v) = (J'(v),u -v) + \b(u -v,u-v), (12) 

where J'(v) = B(v, ■) + 1 is the derivative of J at v. Hence, the first term on the right 
is the first variation of J and it vanishes for a minimizer v = u. Taking u = u^ N ' we 
see that for strictly positive definite B: 



\u 



W -u\\ 2 <- Bin W -u,uW-u) = -( J(u W) - J(uj) y 0. (13) 
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Thus, the Ritz-Lagrange solutions converge to the minimizer even by norm. More 
general conditions for convergence by norm are given in [Til 6. 2 A]. 

Theorem [1] mostly justifies the Ritz-Lagrange method used to solve Problem 1. 
The required properties of J (it) = j$ \{u x ) 2 + fudx and the boundary operator 

rp 

r(u) = (it(0),u(L)) are easily verified, except for the strict convexity. That one 
follows from the Poincare-Friedrichs inequality [T7J 1.6]. We will postpone the discus- 
sion of one other issue until the next section, and focus here on the reason why the 
Ritz method with boundary terms did not work. At first glance, it seems to differ 
from the Ritz-Lagrange method only in the manner the system for «W is derived. 
We would indeed get the same system with the same solutions if we used Eq.Q as 
the starting point rather than Eq.fllJ). The innocent step of replacing Lagrange mul- 
tipliers with their values seems to have spoiled the outcome (the reader is welcome 
to stop reading here and think about this for a moment). 

It was not so innocent after all. It helps to use the hindsight again. We have 
Ao = u x (0) and A^ = —u x (tt) for the exact solution, but for any finite N the trial 
solutions have «£ (0) = u x N \n) = 0, while the Lagrange multipliers are Ao = A^ = 
— ~ 7^ 0. We see that the values of derivatives of the trial solutions at the ends of the 
interval do not converge to the values of the Lagrange multipliers, even though the 
two are equal for the exact solution. Substitution of the limit values in Eq. (J3J) relies 
on exactly the same kind of hidden continuity assumption that was made about the 
residual in the Lanczos tau method. 

How does this reconcile with convergence of the Ritz-Lagrange trial solutions to 
the exact solution? We do have convergence but even for strictly positive definite 
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quadratic functionals it is not strong enough to provide such continuity. Indeed, 
convergence by norm in W\ means that derivatives u x N ^ converge to u x in L2, not 
pointwise. In fact, we will always have u x (0) = Ux (n) = regardless of the 
variational problem as long as we use cosines for trial functions. Unless it so happens 
that the first derivatives of the exact solution vanish at the endpoints there will never 
be convergence of the end point derivatives to the Lagrange multipliers. Lagrange 
multipliers have to be kept as variables and can not be assumed to satisfy relations 
that hold for the exact solution. 



4 Completeness 

Nothing we discussed so far explains why the Ritz-Lagrange method did not work for 
Problem 2. For this problem U = iy|([0, n]), the Hilbert space of functions square in- 
tegrable with their first and second derivatives, with the norm := (IMIl 2 + ||w x ||^ 2 + \\u. 

The functional is J(u) = §^ \{u xx ) 2 — fudx and the boundary operator is again 
T(u) = (w(0), u(L)) T . The space U consists of functions from Wf([0,7r]) that vanish 
at the ends, this space is denoted W| ([0,7r]) in [T7J II. 6]. Problem 3 uses the same 
functional but T(u) = (u{0),u(L),u x (0),u x (L)) , so U = W 2 2 ([0,7r]), the space of 
W| functions that vanish at the ends along with their derivatives. One can check 
that in both cases J and T satisfy all the conditions of Theorem [TJ And yet for 
Problem 2 we did not obtain the exact solution in the limit. 

After inspecting the theorem closely the reader will notice one more condition that 
we did not address yet, because (perhaps) it seemed to be obviously satisfied. But 
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"when you have eliminated the impossible, whatever remains, however improbable, 
must be the truth". That would be the completeness of trial functions, and [H] 
deserves credit for highlighting how far from obvious it can get. Theorem [1] requires 
W 2 completeness in this case. What is obvious, or at least well-known from the 
standard theorems on the Fourier series, is that any function in L 2 ([0,tt]) can be 
approximated by cosines (sines) in the L 2 norm. This is because any L 2 function on 
[0, 7r] can be extended by evenness (oddness) to [— tt, 7r] while remaining in L 2 . But 
that is not even enough for Problem 1, where W 2 completeness was required (this is 
the one issue we skipped). Fortunately, W 2 completeness of cosines reduces to the 
L 2 completeness of sines. 

Lemma 2. The system cos rax, n > is complete and minimal in W 2 ([0, 71*]). 

The minimality above means that the system becomes incomplete after deleting any 
function, the proof is in the Appendix. 

Not so for W%, cosines are incomplete in W|([0, %]). As observed in [13], the 
second derivatives 0, — cos x, —4 cos 2x, . . . , — n 2 cos nx, ... do not include a constant, 
and therefore can not approximate the second derivative of x 2 in L 2 . But then 
cosines can not approximate x 2 in W 2 since its norm incorporates the L 2 norm for 
second derivatives. In [14] the authors add x 2 to the cosine system, but... we shall 
see that adding x 2 is not enough. Here is a quick way to see it. In contrast to 
W 2 l , pointwise values of first derivatives are well-defined and continuous on W 2 2 (this 
follows from the Sobolev embedding theorems [171 1-8], PS 1-8]). Since cosines satisfy 
M x(0) = u x (ti) = any function that can be approximated by them must satisfy the 
same equalities. These are two independent conditions, so the codimension of cosines' 
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linear span is at least two, while x 2 spans only one extra dimension. It follows from 
the next Lemma that the codimension is exactly two and x also needs to be added. 

Lemma 3. The system x,x 2 ,cosnx, n > is complete and minimal in Wf([0, tt]). 

It may seem odd that we have to add x on top of x 2 , after all there is no second 
derivative issue with it. However, while L 2 approximation of the second derivatives 
is necessary for W 2 approximation of functions themselves, it is not sufficient. We 
need to approximate the second derivative and the function itself using the same 
expression. The second derivative of x can most certainly be approximated by cosines 
in L2, as can be x itself, but not while the former approximation is the second 
derivative of the latter. In particular, the cosine series for x diverges in L 2 after two 
differentiations, see Eq.([7]). Note that verifying completeness in a correct space can 
not be avoided even if one uses the usual Ritz method with trial functions satisfying 
the essential boundary conditions. 

What Problem 2 demonstrates is that solutions can not always be approximated 
in the space where the trial functions are complete. Only completeness in the norm 
dictated by the variational functional counts. Completeness of trial functions in a 
norm weaker than the energy norm does not simply weaken the convergence to the 
exact solution, trial solutions may not converge to it at all. Following Lemma [3J 
we should add x and x 2 as the trial functions and represent the trial solutions as 
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(AT) 



b\x + b 2 x 2 + + a n cosnx. The Lagrange functional becomes 



n=l 



£(«<">) = - ^b 2 - |a + 2nbl + ^ 



n=l 



+ A (y + $>n)+A. (^ + ^ + ^ + ^(-1)^) , (14) 



71=1 



2 



n=l 



and the Ritz-Lagrange system is 



f 2 

-y =0 

7T 3 

-—+4^2+77% = 
o 



7rn 



a„ + A + (-1)^ = 



JV 



71=1 



-7T + A + A^ = 7r6 1 + 7r 2 6 2 ^ + ^(-l) n a n = 



JV 



n=l 



(15) 



2 ^ _|_ /■ ]\r 

The equations with A-s immediately yield A = A,r = |, b 2 — — and a n = — 



r?; J 



n > 1. From the second equation in the second column 



N 1 , I t\n i W2J 4 

„ _„(*) ._ 2 y 1 + (- 1 ) _ 1 y 1 

a _a .-Z^ n 4 - 4 2^ A: 4 ^^360 

71=1 fc=l 



as before. Finally, subtracting the second equation in the second column from the 



third we have b± = — irb 2 = Sr. Thus, 



u 



(00) 



o n 4 oo 
(X) = X X H > 

24 24 720 ^ n 

n=l 



3 2/4 3 

7T 7T X X X 2 > 

cosnx = — x — —/. x H 7T Yi\/— 

1 24 M 24 12 /24 



matches the exact solution as expected. 
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This does it for Problem 2, but now we get an unexpected puzzle in Problem 
3. How come we got the correct answer while using an incomplete system of trial 
functions? If a system {0j} is incomplete in U the arguments leading to Theorem [TJ 
still apply to the closure of their linear span in place of U. Therefore, the Ritz- 
Lagrange approximations converge to the minimizer on := U D U$ rather than 



U. For quadratic functional there is a simple relation between u and u<j>. By Eq. (fl2l) 
we have J{u) = J(u) + —B [u — u, u — u) and minimizing J on tl^ is equivalent to 
minimizing B(u — u, u — u) there. As B is a quadratic form, the solution is well known 
to be the orthogonal projection of u to with respect to the inner product defined 
by B [201 IV. 11]. Recall from Lemma [3] that the only functions missing from the 

4 S 2 

%Xj cc 

system of cosines in W$ are x and x 2 . The exact solution u(x) = ir h tt 2 — 

2 v ; 24 12 24 

is S-orthogonal to both of them, despite the misleading appearance of x 2 in the 

formula, as one can check by direct integration. Thus, the missing functions are 

not needed to approximate u simply because u happens to be in (the W 2 closure 

of) the linear span of cosines. This can also be seen directly from its cosine series 

u = > : cosnx, which converges not only in L 2 but even in W$. 

720 ^ n 4 ' y 2 2 

n=l 



5 The Ritz-Lagrange method in higher dimensions 

The Ritz-Lagrange method described in Section [2] can not be applied to multidimen- 
sional boundary value problems. In this section we will develop a suitable general- 
ization and prove that it works. The main distinction is that the boundary operator 
r : U — > V no longer maps into a finite dimensional space. Indeed, in dimensions 
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two and higher the boundary values are not arrays of numbers, but functions on the 
boundary forming an infinite dimensional space V. The induction proof of the key 
Lemma Q] no longer works and its claim itself is false. It is easy to find complete 
systems of functions with no (finite) non-trivial linear combinations satisfying the 
boundary conditions. If we are committed to using arbitrary complete systems of trial 
functions we must find a way to form their linear combinations that satisfy essential 
boundary conditions "approximately" . 

To this end we will use a complete system {if)j} of linear functionals on V, i.e. 
elements of the dual space V* (as with U, the reader may assume that V is a Hilbert 
space, in which case V* = V). If (tpj,Tu) = for all j then Tu = and u e U, so 
we can think of operators T s (u) = ((ipi, Tii), . . . (ip s , Tu)) T as approximations to T, 
and the corresponding spaces U s := {u E U | T s {u) = 0} as approximations to U. 
Assuming V is continuous T s will be also and we can apply Theorem Q] with T s in 
place of T for each s. This gives us a sequence of approximations ui N ^ converging to 
an exact minimizer u s of J on W s . The remaining question is whether we can count 
on u s to approximate the overall minimizer w of J on W. Before proceeding let us 
describe the approximating procedure that our approach suggests. 
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Multidimensional Ritz-Lagrange method. 

To minimize a functional J : U — > K subject to essential boundary condi- 
tions = with r : U — > V select internal trial functions fa , . . . , (j)^ £lA 
and boundary weight functions ipi , . . . , ip s G V* with iV 3> s. A Ritz- 

jV 

Lagrange trial solution u[ N ^ = aifa is obtained by solving the system 

i=i 

of A" + s equations with A^ + s unknowns a\, . . . , a^, Ai, . . . , A s consisting 
of A" internal equations = and s boundary equations (vf)j, Tu) = 0, 
where £ := J(us ) + (X^ s \ Tui ) is the Lagrange functional and := 
Ai^i + • • • + K^Ps is the Lagrange multiplier. 

For higher order equations it is more natural to use several boundary operators 1^ in- 
stead of a single one, e.g. one for function values, another for normal derivatives, etc. 
One can always formally wrap them into a single operator with combined codomain, 
so no additional theoretical discussion is necessary. 

A justification of our method is based on Theorem [1] and Theorem [2] below. The 
reader not interested in justification may skip the rest of this section and look at 
applications in the next one. One piece of bad news is that we can not expect u s 
to converge to u in the same generality as in Theorem [IJ The root cause is that 
the minimizer in U is approximated by elements outside of U, which is why we 
need J to be well-behaved on the entire U, not just U, something avoidable if fa 
do satisfy the boundary conditions. In particular, the values J(u s ) are potentially 
smaller than J(u) because they are obtained by minimizing J on larger subspaces 
U s D Li. As a consequence, standard properties of convex functionals, that we 
relied on in Theorem (U no longer guarantee convergence of J(u s ) to J(u) if u s 

25 



converges only weakly (in technical terms, continuous convex functionals are weakly 
lower semi- continuous, but usually not weakly upper semi-continuous [HI HI. 8. 5], 
[T3l 41.2]). To make our proof work we need to assume a stronger form of convexity. 
For Gateaux differentiable functionals J convexity is equivalent to monotonicity of 
their derivatives, i.e. (J'(u) — J'(v),u — v) > for all u, v [121 H.5.3]. This is a 
generalization of a familiar fact that convex functions have monotone derivatives. 
We will need a form of uniform monotonicity, cf. [121 V.18.6], namely 

(J'(u)-J'(v),u-v)>c(\\u-v\\), (16) 

where c (t) is a continuous monotone increasing function with c(t) = 0. The point is 
that if c( || -Us ||) — > then ||w s || — > and hence u s — > by norm. We are now ready 
to state the main theorem of this section. 

Theorem 2. Suppose J : U — >■ M. is Gateaux differentiable and T : U — >■ V is a 
bounded linear map. Let {ipj} be a complete system in V* and set U s := {u G 
U | (ipj,Tu) — 0, 1 < j < s}. If J grows at infinity on some U SQ , and its derivative 
on it satisfies j?g.(TT6l) then it has a minimizeru on U, as well as minimizers u s on 
all U s with s > sq, and u s > u. The values of J converge to its minimum on U. 

s— >oo 

In examples it is typical that J does not satisfy Eq. (fT6|) on the entire space U, 
but does satisfy on subspaces much larger than U, such as U s . Let us discuss the case 
of quadratic functionals J{u) = \B{u ) u) + l{u) in more detail. From the calculation 
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in Eq.(fl2l) we know that J'{v) = B(v, ■) + I. Therefore, 

(J'(u) — J'(v),u — v ) = B(u, u — v) — B(v, u — v) = B{u — v, u — v) . 

We need B(u,u) > £\\u\\ 2 with e > to satisfy Eq. (fl6l) . i.e. we need B to be 
strictly positive definite. Take the multidimensional analog of the functional from 
Problem 1 as an example, J{u) = J n |(V-u) 2 + fudx, where Q is a domain with 
smooth boundary. It follows from the Poincare-Friedrichs inequality [T71 1.6], [19] 
that B(u,u) = J n |(Vtt) 2 dx is strictly positive definite on W^ity, but it most cer- 
tainly is not on the entire W^^l) since B(u, u) = for any u = const. Nevertheless, 
it still follows from the calculus of variations that B satisfies Eq. (fT6|) on any sub- 
space complementary to the constants, see e.g. [20, VI. 1]. Similar considerations 
apply to other quadratic forms related to the strongly elliptic equations like the bi- 
harmonic equation. They are usually strictly positive definite on complements to 
finite-dimensional subspaces that they annihilate [T9J 22.11]. 

It is worth stressing that Theorems U\ [H do not imply that the double sequence 
ui N ^ converges to u. In fact, let iV < s and let the functionals r*^i, • • • , r*^, where 
r* : V* — > U* is the dual of T, be linearly independent. Then the boundary equations 
(ifjj, Tu) = u) = alone are enough to force ui N ^ = no matter how large iV 

and s are. In practice, this means that one should always take many more internal 
trial functions than the boundary ones, hence iV 3> s. This way for large iV the 
approximation ui will be close to u s by Theorem [TJ while u s in turn will be close 
to u by Theorem [2] if s itself is large enough. 

As in one-dimensional examples one will have to verify completeness of trial func- 
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tions, both internal and boundary, in the appropriate space. One has to be extra 
careful with functionals involving higher order derivatives because the values of func- 
tion and their derivatives have to be approximated simultaneously. Natural spaces to 
use are W£(fl), the spaces of functions with integrable p-th powers along with all of 
their derivatives up to order k. A generalization of the Weierstrass theorem implies 
that polynomials form a complete system in (Q) for any p > 1, any bounded 
domain Q, and any positive integer k (in fact, polynomials are even uniformly com- 
plete [201 II.4.3] ) . However, polynomials may not always be convenient in a particular 
problem. The following Lemma can be useful in finding other complete systems. 

Lemma 4. Let {<pi} and {<pj} be complete systems in W£ (ft) and (Q) respectively, 
where Q and Q are some bounded domains. Then the system {facfrj} is complete in 

w£(n x h). 

If one starts from one- dimensional systems the lemma will only produce complete 
systems in box- like domains [a\, b\] x • • • x [a m , b m }. However, any system of functions 
complete on a domain will be complete on any of its subdomains, so for an arbitrary 
domain one can always use a system complete on the smallest box that contains 
it. Unfortunately, such product systems can not be expected to be minimal even if 
the original systems were. A more targeted choice is to take eigenfunctions of an 
operator on the same domain that is simpler than the one involved, but is somewhat 
similar to it. Various spectral theorems often ensure completeness of eigenfunctions 
in suitable Sobolev spaces [191 22.11a]. 
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6 Multidimensional examples 

In this section we illustrate the multidimensional Ritz-Lagrange method developed 
in Section [5] by applying it to some typical problems. Since calculations by hand 
quickly become intractable we performed them using a computer algebra system. 

Problem 4. Consider a boundary value problem for the Laplace equation V 2 u = 
f in Q, where Q is the unit disk, with the boundary condition u = on dQ and 
/ = cos(a/£ 2 + y 2 )- This equation describes the transverse deflection of a membrane 
fixed everywhere at the boundary and subjected to pressure given by / [201 IV. 10. 3]. 

The profile of / was chosen so that the problem has an analytic solution which 
is not a polynomial. Specifically, one can represent the exact solution as a rapidly 
convergent series 

u(r) = 7 + oos(l) - Ci(l) - cos(r) + £ r ^ (17) 

where r = x 2 + y 2 , 7 := lim^oo ^ Ylk=i i — mn ) ^ s ^ ne Eiiler-Mascheroni constant, 
and Ci(x) := — J°° dt is the cosine integral. 

To solve this problem we use the multidimensional Ritz-Lagrange method. The 
variational problem is to minimize the functional J{u) — f - (Vu) 2 + fu dxdy, which 
gives the total potential energy of the membrane, subject to the boundary condition. 
In the notation of Section [5] we take U = Wl(VL) with V being the restriction of u to 
the boundary dQ. Moreover, T is continuous if we take V = L 2 (dVt). Our internal 
trial functions are the monomials, which obviously do not satisfy the boundary con- 
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N N 

dition, and the trial solution is u^ N ^ = x t ~ 1 j/~ 1 . Note that A" of Section 

i=\ 3=1 

[5] will be A^ 2 here because of double indexing. As the boundary weight functions 
we choose the piecewise linear ones on uniform partitions of dVt. Unlike the usual 
choices for a circle, e.g. the trigonometric functions, these ones can be used on a wide 
variety of boundaries. Instead of using a single indexed system ip^ it is convenient to 
split it into the constant g% and the linear hk parts. If the boundary is partitioned 
into s segments we have 



,„.(«): { « s and h t (0):=l " s 

0, otherwise 0, otherwise. 

(18) 

Therefore, the number of boundary weight functions, denoted s in Section [5j will be 

s-l 

2s here. The Lagrange multiplier has the form = (ckgk(d) + dk hk(9)), and 



fc=0 



the Lagrange functional is £(u^ J = J(u^) + / X^u^ N ' da. The unknown coef- 

^ ' Jan 



and 2s boundary / gkU^ da = h^u^ da = equations. 
Jan Jan 



ficients a, ,-, Ck and dk are now determined from the system of A^ 2 internal — — = 

da 



l :J 



'an Jan 

The relative errors of the Ritz-Lagrange solutions versus the exact solution u 



(JTTJ) are shown in Table [T] as the percentages of the maximum deflection at x = 
and y = . They are quite small considering that one has to determine A^ 2 + 2s 
coefficients in each case. Note that we always keep A^ 2 > 2s as recommended in the 
description of the method to ensure that the system matrices have full rank and are 
invertible. 
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Central Error % 


Boundary Error % 


N = 3,8 = 2 


-4.04 


-1.33 


N = 4, a = 3 


-4.05 


-0.1 


A^ = 5, s = 4 


-0.04 


0.03 



Table 1: Central (x — 0, y — 0) and boundary (x = 0,y = 1) error relative to the 
maximum deflection of the membrane of unit radius. 

Problem 5. Consider a problem of bending a uniformly loaded, simply supported 
on all sides (SS-SS-SS-SS), isotropic, square plate of constant thickness, unit stiffness 
and unit edge length. Simply supported means that u = on dQ. We do not need 
to list the natural boundary conditions since a variational formulation incorporates 
them automatically. The variational functional giving the potential energy of the 
plate is [20, IV.10.3]: 

f 1 f 1 1 (fd 2 u\ 2 fd 2 u\ 2 d 2 ud 2 u , ,fd 2 u\\ ril 

J(u) = /o I + W) +2v ^W + 2il - v) {^y))- fudxdy 

(: 

where u is the displacement of the plate, v is the Poisson ratio, and / is a distributed 
load. 

The Euler-Lagrange equation induced by Eq. ffl9|) is the biharmonic equation 
V 2 V 2 m = /, the terms multiplied by the Poisson ratio form a divergence and only 
affect the natural boundary conditions. As the internal trial functions we choose the 
products of cosines X^x) = cos((i — 1)ttx) and Yi(y) = cos((z — l)7r?/), so that the 
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trial solution is = ^] Cjj Xj-\ (x) Yj-\ (y) ■ Obviously, the trial functions do 

i=l j=l 

not satisfy the boundary condition. The Lagrange functional is 



Z{u^ = j(u^ + j\[ s \x)u(x,0)dx 

+ [ \¥{y)u{l,y)dy- f X^\x)u(x, 1) dx - [ \[ s \y)u(0,y) dy , (20) 
Jo Jo Jo 

where for convenience we split the Lagrange multiplier A^ into its restrictions A- 5 '* 
to each edge of the plate. This way we can represent the set of the boundary 
weight functions as the union of four sets selected separately for each edge, namely 

s 

\[ s \z) = y^Aijeos ((j — 1)ttz), where Ajj are the unknown coefficients. The num- 

ber of internal equations here is again N 2 , and the number of the boundary equations 

is 4s, so the non-degeneracy condition is iV 2 > 4s. 

The exact solution to this problem can be expressed as a rapidly convergent series, 

we use its first ten terms to calculate the errors. We do not tabulate them here, 

because they are very large (up to 70%), and increasing the number of terms does 

not improve the approximation. At this point, the reader should not be surprised, 

indeed we are dealing with the same mistake as in the one-dimensional Problem 2. 

As we pointed out in Section EJ the system of cosines is incomplete on the interval, so 

the system of their products naturally is incomplete on the product of intervals that 

represents the plate. A more down to Earth explanation is that products of cosines 

Ou 

have vanishing normal derivatives — — on all edges and all their linear combinations 

on 

inherit this property. This does not matter for second order equations because the 
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normal derivatives are discontinuous on the relevant spaces, but it does matter for 

the higher order equations like the biharmonic equation. 

One can check that in the weak formulation of Problem 5 because of the vanishing 

normal derivatives the boundary terms that multiply the variation of the solution's 

derivatives get removed, so we ended up solving a different variational problem. 

Indeed, the choice of cosine products unwittingly enforces an additional boundary 
Ou 

condition, — = on d£l. Together with u = on <9f2 this describes, physically, a 
on 

plate clamped on all sides (C-C-C-C) rather than a simply supported one. Thus, we 
should be comparing our Ritz-Lagrange solutions to the answers for the C-C-C-C 
plate (cf. Problem 3 in the one-dimensional case). Unfortunately, an analytic solution 
for a plate clamped on all sides is not known, so we used the values obtained in [2U 
VI. 44] to make the comparison. The relative errors as percentages of the maximum 
deflection at the center of the plate are shown in Table [2J 





Central Error % 


Boundary Error % 


N = 4,s = 2 


52.76 


-50.92 


N = 6,s = 3 


-3.3 


-1.22 


N= 8,3 = 4: 


0.063 


-1.39 


N= 10, s = 5 


-0.12 


-0.09 



Table 2: Central (x = 0.5, y = 0.5) and boundary (x = 0, y = 0) error relative to 
maximum bending deflection of an C-C-C-C square plate of unit size. 

Section [3] also gives us a way to solve the original problem, we just need to 
complete the system of cosine products. By Lemma H] it suffices to complete the 
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cosines on the interval and take the products from the completed system. Namely, 
we take the products of x, x 2 , cos((z — l)irx) and y, y 2 , cos((i — l)vn/) as the new trial 
functions, and keep the rest of the above setup intact. The relative errors for the 
Ritz-Lagrange solutions with the completed system against the known series solution 
[21 8.2.4] are shown in Table [3] and demonstrate the validity of the method. 





Central Error % 


Boundary Error % 


N = 5,s = 2 


26.14 


-89.03 


N = 6,s = 3 


2.81 


-3.7 


N= 8,3 = 4: 


1.68 


-4.4 


N= 10, s = 5 


0.56 


-1.1 



Table 3: Central (x = 0.5, y = 0.5) and boundary (x = 0,y = 0) error relative to 
maximum bending deflection of an SS-SS-SS-SS square plate of unit size. 

As a final demonstration, we apply the Ritz-Lagrange method to a boundary 
eigenvalue problem for square plates. The eigenmodes describe standing vibrations 
of a plate, and their zeros (nodal curves) are known as Chladni figures [221 5.1]. The 
problem has attracted a lot of attention from both analytic and numerical viewpoints, 
indeed Ritz himself applied his method to it in his original paper. Boundary eigen- 
value problems are somewhat beyond the scope of the theory in Section [51 which 
deals with linear constraints only, because under the Rayleigh-Ritz approach one 
needs to impose an additional quadratic normalization constraint | f n u 2 dx = 1 to 
solve for the eigenmodes [191 18.5], [20l VI. 1.1], [221 5.2]. However, the general ap- 
proach of Section [5] remains valid, and with some extra effort one can justify applying 
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the Ritz-Lagrange method to problems with non-linear constraints along the same 
lines. 

Problem 6. Consider a uniformly loaded, simply supported on all sides, isotropic 
square plate of constant thickness with unit edge length. The potential energy J of 
the plate is given by Eq. (ll9j) without the distributed load term. The boundary 
eigenvalue problem can be interpreted as finding extrema of J(u) subject to the 
boundary condition u = on dQ, and the normalization constraint \ J n u 2 dx = 1. 

Compared to Problem 5 the Lagrange functional acquires an additional term 
fi(f n u 2 dx — 1) and an additional equation, which amounts to the normalization 
constraint on the eigenvectors. Of course, in practice one can ignore this equation 
and simply use standard methods for finding eigenvectors. We keep the choices for 
the internal and the boundary weight functions from Problem 5. Let c denote the 
vector of internal coefficients and A denote the vector of boundary coefficients 
Xij. In terms of c and A the Lagrange functional can be conveniently represented 
as = \cFKc — ii ^c T Mc — l) + \ T Lc, where K and M are matrices of size 

TV 2 x N 2 obtained by integrating the internal trial functions, see 8.2.7], and L is 
a 4s x N 2 matrix obtained by integrating the boundary weight functions. Matrix 
L can be obtained by multiplying the boundary equations with the corresponding 
Lagrange multiplier functions and extracting the coefficients of Cjj and \ t j after the 
integration. We note that the boundary equations can be written as Lc = 0. Finally, 
differentiating the Lagrangian with respect to Cy and \j we are led to the following 
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generalized eigenvalue problem: 



K 


L T 


L 






M 












\ / \ 

c 



/ 



(21) 



For this eigenvalue problem to be solvable one needs L to have the maximal rank 4s, 
which is ensured by the non-degeneracy condition N 2 3> 4s as in Problem 5. The 
eigenvalues \ii = uf approximate the squares of the natural frequencies of the plate's 
vibrations. 

With N = 10 and s = 5 we obtain a set of approximate non-dimensional natural 
frequencies Ui, first nine of which are shown below. Since the eigenmodes are known 
to be of the form sin(7rma;) sin (7m?/) we change the single index notation to u mn and 
arrange the frequencies in a square pattern 



19.61797 49.06479 98.33527 
49.06479 77.55724 126.62091 
99.03778 126.62091 177.73211 



(22) 



The exact values are taken from [21 8.2.4], repeated frequencies corrrespond to mul- 
tiple eigenvalues with the eigenmodes symmetric along different axes: 



19.73920881 49.34802202 98.69604404 
49.34802202 78.95683523 128.3048573 
98.69604404 128.3048573 177.6528793 



(23) 
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One can see that the estimated frequencies are slightly lower than the exact ones. 
This is in contrast with the application of the usual Ritz method, where the estimated 
frequencies are always higher. From a physical viewpoint, the latter happens because 
replacing an infinite system with a finite one is equivalent to imposing additional 
constraints, which tend to raise the stiffness of the system, and hence the frequencies. 
This assumes however that all the boundary constraints are enforced in both systems, 
i.e. that the trial functions satisfy the essential boundary conditions. 

In the Ritz-Lagrange method the trial functions do not satisfy the essential con- 
ditions, and even the trial solutions are forced to satisfy them only approximately. 
In other words, i.e. we are effectively relaxing the boundary constraints in addition to 
imposing additional ones through discretization. This relaxation lowers the frequen- 
cies (because a plate with fewer constraints is less stiff) and counteracts the effects of 
discretization. If we were able to impose the boundary conditions everywhere along 
the boundary the estimated frequencies would have been higher than the exact ones 
just as in the usual Ritz method. 

From a mathematical viewpoint, this effect is also natural since the eigenvalues 
are the minima of a quadratic functional on subspaces of the original space [201 
VI. 1.1]. In the Ritz-Lagrange method we approximate them by using functions from 
a larger space (by relaxing the boundary conditions), thus lowering the minima that 
can be attained. In particular, one can see from the proof of Theorem [2] that the 
values of the functional at the approximating elements are potentially smaller than 
at the sought minimizer. 
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7 Conclusions 



We developed a general extension of the Ritz method to systems of trial functions that 
do not satisfy the essential boundary conditions, and proved its convergence. The 
method is based on treating the essential conditions as variational constraints and 
removing them using the Lagrange multipliers. Here are some general observations 
of the workings of the method. 

• The variational functional has to be well-behaved not only on the energy space 
of the problem, but on its extension that contains the trial functions. Suf- 
ficiently good behavior is a strong form of convexity, which in the case of 
quadratic functionals means that the boundary value problem is strongly ellip- 
tic. 

• The systems of trial functions must be complete in the norms consistent with 
the functional, which usually restrict to the energy norms on the energy space 
of the problem. Although similar requirement applies to the usual Ritz method, 
here it is much easier to encounter systems that appear complete but are not 
due to effects at the boundary. 

• The Lagrange multipliers have to be treated as additional variables in the 
approximating systems. They can not be eliminated by substituting the trial 
solutions into the variational formulas for them in terms of the exact solution. 
These formulas are discontinuous in the relevant norms. 

• In multidimensional problems the boundary conditions incorporate infinitely 
many constraints, and to obtain a finite dimensional approximating system one 
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has to select boundary weight functions in addition to the trial functions. The 
number of trial functions has to be significantly larger than the number of the 
boundary weights, otherwise the approximating system may be inconsistent or 
only have the trivial solution. 

• In multidimensional problems the approximating values of the functional may 
approach the exact value from below rather than above, as in the usual Ritz 
method, because the minimization takes place on a larger space of functions 
not satisfying the boundary conditions exactly. 

• The method can be applied to boundary eigenvalue problems interpreted along 
the lines of Rayleigh-Ritz as minimization problems on subspaces of the orig- 
inal space with the additional normalization constraint. Due to the presence 
of the Lagrange multiplier variables the resulting finite dimensional problem 
is a generalized eigenvalue problem (A — fj,B)x = instead of the ordinary 
one with B = I. In multidimensional vibrational problems the approximate 
eigenfrequencies obtained in this way may be lower than the exact ones, in 
contrast to the Ritz method where they are always higher, due to relaxation 
of the boundary constraints. 

As is well-known pQE], the Ritz method leads to the same approximating systems 
as the Galerkin method, but the latter can also be applied to non-optimization prob- 
lems. It is quite intriguing whether the Ritz-Lagrange method developed here can be 
extended to a 'Galerkin-Lagrange method' for non- variational problems. There is no 
Lagrange functional to be had in such problems, but one can formally add boundary 
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terms multiplied by extra variables ('Lagrange multipliers') to a weighted residual of 
the problem. An approximation of the boundary conditions should also be added to 
the usual Galerkin system. However, as we saw in the Lanczos tau example such a 
straightforward approach is not likely to work. Indeed, in the Ritz-Lagrange method 
we add the Lagrange boundary terms not to the bare weighted residual, but to an 
integrated by parts expression with some boundary terms of its own. This suggests 
that in a correct generalization the problem has to be rewritten in a weak form J2j 
7.5.1] before the Lagrange-like boundary terms are added. 

Another complication is the role of the natural boundary conditions. In varia- 
tional problems they are enforced automatically, so there is no point in adding them 
as constraints and introducing additional Lagrange multipliers. An example in [H] 
even shows that attempting to do so leads to worsening the convergence of the trial 
solutions. However, it is unclear how the natural conditions should be enforced in 
a 'Galerkin-Lagrange method'. Still, the distinction between the natural and the 
essential boundary conditions in [31 1.1.2] (by the order of the derivatives entering 
them) makes sense even for non-optimization problems, and it has been shown in p[] 
that in some cases the weak form of a problem provides an enforcement mechanism 
for the natural conditions. 

Appendix: Proofs 

Proof of Lemma Q] . Let be a complete system in U. Since V has an s-dimensional 
image we can represent it as Yu = {(Ti, u), . . . , (r s , u)) T , where Tj are bounded 
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linear functionals. Assume without loss of generality that they are linearly inde- 
pendent, otherwise some of them can be dropped without changing U. Set Uk : = 
{ueU \ (Ti, u), ... , (r^, u) = 0}, we will construct a complete system in each Uk by 
induction on k. Since U = U S the process concludes in s steps. 

For k = 1 we must produce a complete system of linear combinations in U\ : = 
{u G U | (ri, u) — 0}. Without loss of generality, (ri, 0i) 7^ since {0j} is complete 

and Ti can not vanish on all 0j. We claim that & := & — ' \ (pi G Wi form the 

(ri,0i) 

desired system. Let u G U\ C W and be the coefficients such that X^ili a i — 
e for a given e > 0. By definition of (p i: 

N N /r v~»JV , v 

= 2^ <Pi (tY^i) — ' 

To estimate the second term we find, 

N N N 

|(ri,y^a^)| = |(ri,y^ai 4>j -u) + (ri,u)| < ||ri|| y^a, <^|| < ||ri||e. 

i=l i=l i=l 

N _ up || ||^ || 

Therefore, \\u — > ajSj 

II < 1 + ,-n , ; e, and since w, e are arbitrary com- 
^ v l( r i,<Ml y 

pleteness of 0j follows. 

Let be a complete system in Uk from the preceeding step. Linear indepen- 
dence of Tj guarantees that r^+i does not vanish on some 0j, which we may as well 
take to be 4>i- Apply the process above with T 1 replaced by T k +i and U\ replaced by 
Uk+i to obtain <f> i . Then <p i are linear combinations of <pi (and hence of the original 
(pi ), belong to Uk+i and are complete in it by the same argument. This concludes 
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U — >OD 



the induction step. □ 

Proof of Theorem [TJ A standard argument from convex analysis shows that if J(u) > 

oo on hi then J has minimizers on hi and there is a weakly convergent subsequence 

u {N k ) u (oo) 

|12[ III. 10.3], [111 6.2]. For convex functionals strong continuity 

k— >oo 

implies weak continuity so J 

( u (N k )\ y J( u (°°)), and moreover u (oo > E hi since 

k— >oo 

= T(u^ Nk A > T(u(°°A. But for large enough iV there isaueW^' arbitrarily 

close to a minimizer w of J on W by Lemma (TJ hence by continuity J{u) is arbitrarily 
close to the minimal value J m i n - But j(u^ Nk ') can not exceed J(u) for > N since 
it is a minimizer on hi^ Nk \ so J m ; n < J(u( Nk ^) < J(u) = J m i n + e. After passing to 
limit we have that J(m < -°°- ) ) = J mm , i.e. -u^ 00 ) is a minimizer of J on hi. 

If we assume additionally that J is strictly convex on hi then m is unique and 
the entire sequence u ( N > (which is now also uniquely defined) converges to it at least 
weakly [11, 6. 2 A]. □ 

The next two proofs use equivalent norms (inner products) on ([0,7r]) and 
WfClP, ft"]) respectively. Two norms are equivalent if they define the same notion of 
convergence, for equivalent norms on Sobolev spaces see [171 1-8] and especially [TBI 
1.9]. 

Proof of Lemma [2] The following inner product is equivalent to the usual one on 

([0, 7r]): (u, v)o :— u(0)v(0) + / u x v x dx. To prove completeness it suffices to 

Jo 

show that any function w orthogonal to all cosines must be 0. For such w we have 

(w, l)o = w(0) = and hence (w, cosnx)o — I w x ■ (— nsinnx) dx = for n > 1. 

Jo 

Thus, w x is L2 orthogonal to sinnx for all n > 1. Since the latter form an orthogonal 
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basis in L 2 ([0, 7r]) we must have w x — a.e. But then by the Fundamental Theorem 

of Calculus w{x) = w(0) + w t dt = a.e. establishing completeness. Being an 

Jo 

orthogonal basis in L 2 ([0,7r]) cosines must be minimal there, and therefore in any 
space with a stronger norm, which includes W 2 1 ([0, 7r]). □ 

Proof of Lemma OS An equivalent inner product on VF 2 2 ([0,7r]) is 

/'TT 

(u, v)o := u(0)v(0) + u x (0)v x (0) + / u xx v xx dx. Consider w orthogonal to all cosines, 

Jo 

PIT 

then we have (w, l)o = w(0) = and (w, cosnx) = / w xx ■ (— n 2 cosnx) dx = for 

Jo 

n > 1 because all sines vanish at 0. In particular, w xx is L 2 orthogonal to cosnx for 
all n > 1. But orthogonal complement of the latter in L 2 consists of constants, so 
w xx = const and w(x) = ax 2 + bx + c. Since w(0) = free term is and w is a linear 
combination of x and x 2 . Thus, orthogonal complement to cosines is spanned by x 
and x 2 proving completeness. 

For minimality notice that by direct calculation (x,cosnx)o = {x 2 ,cosnx)o = 
(x,x 2 )o = 0, i.e. x and x 2 are orthogonal to all cosines and to each other. This 
means that neither one of them can be deleted without loosing completeness. It also 
means that if a cosine can be approximated in W2 by other cosines combined with 
x and x 2 then it can already be approximated by other cosines alone. But the latter 
can not be done with arbitrary precision even in L 2 , let alone in W 2 . □ 

Proof of Theorem [2j Since Li C • ■ ■ C U2 C Ui and the minimum on a larger 
space can not get bigger we have J(u) > • • ■ J(u 2 ) > J(ui)- Thus, the numerical 
sequence J(u s ) is bounded. Moreover, u s e U so for s > s , so ||w s || < M < 00 
for s > So since J grows at infinity on U So . Recall that u, u s &re the minimizers 
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of J on U, IA S respectively, and therefore the derivatives J'(u), J'(u s ) vanish when 
paired with elements of the corresponding subspaces. In particular, (J'(u),u) = 
(J'(u s ),u) = (J'(u s ),u s ) = and 



(J'(u s ) - J'(u),u s - u) = -(J'(u),u s ). 



(A.l) 



We will prove that the last expression converges to when s — > oo. Condition Eq. ffTo 7 ]) 
then implies c (u s — u) > and hence u s > u as claimed. Convergence of J(u s ) 

s—>oo s—>oo 

follows from continuity of J. 

We now prove convergence in Eq.f lA.ll) . Since u is a minimizer on hi the functional 
J'(u) vanishes on any element from it. The subspace of functionals that vanish on 
the entire U = {u E U | Tu = 0} is the closed linear span of {T*ipj} in U* . Indeed, 
if {r*^} did not span it there would exist, by the Khan-Banach theorem, a u such 
that (r*ipj,u) = (ipj, Tu) = for all j, while Tu ^ 0, contradicting the completeness 
of {ipj}- Thus, for any e > there exists a linear combination £ = Y^j=\ a j^*' l l J j suc h 
that || J' (u) — £|| < e. But then £ e U n , and for s > n we have (£,u s ) = 0, so 

\(J'(u),u s )\<\\J'(u)-a \\u s \\<Me. 
Since e is arbitrary (J'(u),u s ) > 0. □ 

s— >oo 

Proof of Lemma [41 In the multiindex notation an equivalent norm on Wp{V) is 
given by 



-EE 

i=0 \a\=i 



E 

\a\<k 



where the L p norm is just ||/|| Lp := (J v \f\ p dQ*. life W£(tt) and / G W*{VL) then 
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it follows from the Fubini theorem that ||//||l p = ||/||l p ||/||l p since / and / depend 
on different variables. Let x and x denote the variables on Q and Q respectively, so 
that £ = (x, x) is the variable on Q x Q. Then we estimate 



\FF\ 



s E 

\0\<k,\~f\<k 



E| 

\a\<k 

\d^F 



d\°\FF 



E 

m+h\<k 



\d\®F 




Q\l\p 


\ dx 13 


Lp 


dx"< 



dx 13 



Q\i\p 



dx"f 



E 

\P\<k 



d\P\F\\ 
dxP IL„ 



E 

|7|<fc 



dx"i 



\\\\ 



k\\F\ 



Since <pi are complete any monomial x 13 can be approximated to any precision e > 
in Wp by their linear combination = ^E a, 0j, and analogously can be ap- 



proximated by a linear combination 



But 



is a 



linear combination of 
arbitrarily small: 



j, while the difference between the products can be made 



x^x 1 — 4>4>\\w$ — ||x (5^ — 4>) + (x — <j>)<ft\\w£ 

< Wx^WwkWx 1 - (f)\\ W k + \\af - </>\\ W k\\$\\ W k <e(\\x^\\ W k + ||x 7 ||^ +e) , 



where the first inequality follows from the above estimate. Hence any product of 
monomials, and therefore any polynomial, can be approximated in by linear 
combinations of <pi<pj. By the generalized Weierstrass theorem J20J H.4.3], polynomi- 
als are complete in Wt (ft x Q), and hence so is the system {(j)^}. □ 
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